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Abstract 

The theory of random matrices with eigenvalues distributed in the complex 
plane and more general "/3-ensembles" (logarithmic gases in 2D) is reviewed. The 
distribution and correlations of the eigenvalues are investigated in the large N 
limit. It is shown that in this limit the model is mathematically equivalent to a 
class of diffusion-controlled growth models for viscous flows in the Hele-Shaw cell 
and other growth processes of Laplacian type. The analytical methods used involve 
the technique of boundary value problems in two dimensions and elements of the 
potential theory. 
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1 Introduction 



Applications of random matrices in physics (and mathematics) are known to range from 
energy levels statistics in nuclei to number theory and from quantum chaos to string 
theory. Most extensively employed and best-understood are ensembles of hermitian or 
unitary matrices (see e.g. [T]-[I]). Their eigenvalues are confined either to the real axis or 
to the unit circle. In this paper we consider more general classes of random matrices, with 
no a priori restrictions to their eigenvalues. Such models are as yet less well understood 
but they are equally interesting and meaningful from both mathematical and physical 
points of view. (A list of the relevant physical problems and corresponding references 
can be found in, e.g., [5].) 

The progenitor of ensembles of matrices with general complex eigenvalues is the sta- 
tistical model of complex matrices with the Gaussian weight introduced by Ginibre [6] 
in 1965. The partition function of this model is 



[D$]exp (- jtr$ f $) . (LI; 



Here [-D$] = Ylij d(dt <&ij)d(5s 3>y) is the standard volume element in the space of N x N 
matrices with complex entries $ij and t is a (real positive) parameter. Along with the 
Ginibre ensemble and its generalizations we also consider ensembles of normal matrices 
[7j, i.e., such that $ commutes with its hermitian conjugate 

Since one is primarily interested in statistics of eigenvalues, it is natural to express 
the probability density in terms of complex eigenvalues Zj = Xj + iyj of the matrix $. It 
appears that the volume element can be represented as 

[D^ ( xl[\z i -z j \ 2 l[d 2 z i . (1.2) 

i<j i 

If the statistical weight depends on the eigenvalues only, as it is usually assumed, the other 
parameters of the matrix (often referred to as "angular variables" ) are irrelevant and can 
be integrated out giving an overall normalization factor. In this case the original matrix 
problem reduces to statistical mechanics of A" particles with complex coordinates z~ in the 
plane. Specifically, the factor Ili<j R ~ z j\ > being equal to the exponentiated Coulomb 
energy in two dimensions, means an effective "repelling" of eigenvalues. This remark leads 
to the Dyson logarithmic gas interpretation [8], which treats the matrix ensemble as a 
two-dimensional "plasma" of eigenvalues in a background field and prompts to introduce 
more general "/3-ensembles" with the statistical weight proportional to Ili<j \ z % ~ z j\ ■ 

It is also natural to consider matrix ensembles with statistical weights of a general 

form, exp ^ — — tr W ($>\ , with a background potential W . An important observation 

made in [9] and developed in subsequent works [TOj, [TTJ [12j [13] is that evolution of an 
averaged spectrum of such matrices as a function of t, as A" — > oo, serves as a simulation 
of Laplacian growth of water droplets in the Hele-Shaw cell (for different physical and 
mathematical aspects of the latter see [TU [15], IE]). To be more precise, in the limit 
A^ — > oo, under some weak assumptions about the statistical weight, the eigenvalues are 
confined, with probability 1, to a compact domain in the complex plane. One can ask 
how its shape depends on the parameter t. The answer is: this dependence is exactly 
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the Laplacian growth of the domain with zero surface tension. Namely, the edge of the 
support moves along gradient of a scalar harmonic field in its exterior, with the velocity 
being proportional to the absolute value of the gradient. The general solution can be 
expressed in term of the exterior Dirichlet boundary value problem. 

This fact allows one to treat the model of normal or complex random matrices as a 
growth problem [12]. The advantage of this viewpoint is two-fold. First, the hydrody- 
namic interpretation makes some of the large N matrix model results more illuminating 
and intuitively accessible. Second and most important, the matrix model perspective 
may help to suggest new approaches to the long-standing growth problems. In this re- 
spect, of special interest is the identification of finite time singularities in some exact 
solutions to the Hele-Shaw flows with critical points of the normal and complex matrix 
models. 

2 Random matrices with complex eigenvalues 

We consider square random matrices <£> of size N with complex entries The proba- 
bility density is assumed to be of the form P(<&) oc e^ trW ^\ where the function 
(often called the potential of the matrix model) is a matrix-valued function of $ and $t 
such that = and H is a parameter introduced to stress a quasiclassical 

nature of the large N limit. The partition function is defined as an integral 

Z N = J[D$]e* tTW W (2.1) 

over the matrices with the integration measure [D<&] to be specified below in this section. 

We consider two ensembles of random matrices $ with complex eigenvalues: ensemble 
C of general complex matrices (with no restrictions on the entries except for det $ ^ 0) 
and ensemble Af of normal matrices (such that [<&, <$>'] — o). 

2.1 Integration measure 

The integration measure has the most simple form for the ensemble of general complex 
matrices: 

N 

This measure is additively invariant and multiplicatively covariant, i.e. for any fixed (non- 
degenerate) matrix A 6 C we have the properties [D(<& + A)] = [D<$>] and = 
[D(A$)] = | det A\ 2N [D$\. It is also clear that the measure is invariant under transfor- 
mations of the form $ — ► WOfU with a unitary matrix U. 

The measure for Af is induced by the standard flat metric in C, \ | 2 = tr (<5 < M$^) = 
l^iil 2 via the embedding Af C C. Here Af is regarded as a hypersurface in C defined 
by the quadratic relations = As usual in the theory of random matrices, one 

would like to integrate out the "angular" variables and to express the integration measure 
through the eigenvalues only. 
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The measure for Af through eigenvalues [TJ [7] . We derive the explicit represen- 
tation of the measure in terms of eigenvalues in three steps: 

1. Introduce coordinates in Af C C. 

2. Compute the inherited metric on Af in these coordinates: ||5$|| 2 = g a pd^ a d^. 

3. Compute the volume element [D$>] = ^\ det g a p \ l\ a d^ a . 

Step 1: Coordinates in Af. For any matrix $, the matrices H% = |($ + H 2 = 
— $t) are Hermitian. The condition [$, $ T ] = is equivalent to [ifi, #2] = 0. Thus 
#1,2 can be simultaneously diagonalized by a unitary matrix U: 

Hi = UXW, X = diag{xi, . . . , xjv} 
H 2 = UYW, Y = dmg{ yi ,...,y N }. 

Introduce the diagonal matrices Z = X + iY, Z = X — iY with diagonal elements 
Zj = Xj + iyj and Zj = Xj —iyj respectively Note that Zj are eigenvalues of $. Therefore, 
any $ e Af can be represented as $ = UZW, where U is a unitary matrix and Z is the 
diagonal matrix with eigenvalues of $ on the diagonal. The matrix U is defined up to 
multiplication by a diagonal unitary matrix from the right: U — > UU^ ag . The dimension 
of Af is thus 

dim (Af) = dim (U) - dim (W diag ) + dim (C diag ) = N 2 - N + 2N = N 2 + N 
(here U C C is the submanifold of unitary matrices). 

Step 2: The induced metric. Since $ = f/Zf/ 1 ", the variation is 5$ = U(5u ■ Z + 5Z + Z ■ 
5iv)W, where 5u^ = U5U^ = —5u^. Therefore, 

||5$|| 2 = tr(<5<M$ f ) = ti (5 Z 5 Z) + 2 ti(5uZ5uZ-(5u) 2 ZZ) 

N N 
j=l j<k 

(Note that Sujj do not enter.) 

Step 3. The volume element. We see that the metric g a p is diagonal in the coordinates 
$1(5 Zj), $s(5zj), $t(5ujk), Q(5ujk) with 1 < j < k < N, so the determinant of the diagonal 
matrix g a p is easily calculated to be | det g a p\ = 2 N ~ N Y\f <k \zi — z^. Therefore, 

N 

[£>$] cx \DU}' \A N ( Zl , . . . , z N )\ 2 J] d 2 Zj , (2.2) 

i=i 

where d 2 z = dxdy is the flat measure in the complex plane, [DU]' = [DU]/[DUdia,g\ is 
the invariant measure on U/U^^ and A N is the Vandermonde determinant: 

N 

A N (z u ...,z N ) = U(zj-z k ) = det N (z^). (2.3) 

j>k 
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The measure for C through eigenvalues. A complex matrix $ with eigenvalues 
Zi,...,z N can be decomposed as $ = U(Z + R)W, where Z = diag{zi, . . . , z N } is 
diagonal, U is unitary, and R is strictly upper triangular, i.e., Rij — if i > j. These 
matrices are defined up to a "gauge transformation": U UU^ &g: R — > R C/diag- It 
is not so easy to see that the measure factorizes. This requires some work, of which the 
key step is a specific ordering of the independent variables. The final result is: 

[DQ] oc [DU]' \Y[d 2 R kl ) \A N (z t )\ 2 {[ d 2 Zj . (2.4) 

\k<l J j=l 

The details can be found in the Mehta book pQ. 



2.2 Potentials 

For the ensemble M the "angular variables" (parameters of the unitary matrix U) always 
decouple after taking the trace trW(&) = J2j W(zj), so the potential W can be a function 
of $, $t of a general form W($) = a nm^ n {^) m - The partition function reads 

r N 

Z N = / lA^z^Ue^^d 2 ^ (2.5) 

J 3=1 

where we ignore a possible A^-dependent normalization factor. From now on this formula 
is taken as the definition of the partition function. 

The choice of the potential for the ensemble C is more restricted. For a general 
potential, the matrix U in $ = U(Z +R)W still decouples but R does not. An important 
class of potentials when R decouples nevertheless is VT($) = — + V(&) + V(&), where 
V(z) is an analytic function of z in some domain containing the origin and V(z) = V(z). 
In terms of the eigenvalues, 

W(z) = -\z\ 2 + V{z) + V(z). (2.6) 

In what follows, we call such a potential quasiharmonic. In this case, tr = tr (ZZ) + 

tr (RR^), tr ($ n ) = tr (Z + R) n = tr Z n , and so 

[ [D$]e^ TW W = Cn [ \A(z i )\ 2 ]Je^ w ^d 2 z k , (2.7) 
Jc J k 

where is an A-dependent normalization factor proportional to the gaussian integral 
j[DR]e-rM^). 

As an example, let us consider the quadratic potential: W(z) = —\z\ 2 + 2^t (tiz+t 2 z 2 ). 
The ensemble C with this potential is known as the Ginibre-Girko ensemble [6j [T7J. In 
this case the partition function (12.51) can be calculated exactly [IS] : 

7 _ 7(0)/, n. |2 A -JV 2 /2__ ( N tfh+t 2 ^ + |ti| 2 ^ 

Z N -Z N (1-A\t 2 \) exp^- j, (2.8) 

N 

where zffl = f^ N +n ^ 2 tt n Y\ k\ is the partition function of the model with W = — \z\ 2 . 

k=l 
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Note that for quasiharmonic potentials the integral (I2.5P diverges unless V is quadratic 
or logarithmic with suitable coefficients. The simplest way to give sense to the integral 
when it diverges at infinity is to introduce a cut-off, i.e., to integrate over a big but finite 
disk of radius Rq centered at the origin. Just for technical simplicity we assume that a) 
V is a holomorphic function everywhere inside this disk, b) W has a maximum at the 
origin with W(0) = and no other critical critical points inside the disk, c) At \z\ = Rq 
the potential W is bounded from above by a constant B < 0. The large N expansion is 
then well-defined. For details and rigorous proofs see [T9 l 1201 I2T] . 

2.3 The Dyson gas picture 

The statistics of eigenvalues appears to be mathematically equivalent to some impor- 
tant models of classical statistical mechanics, with the eigenvalues being represented as 
charged particles in the plane interacting via 2D Coulomb (logarithmic) potential. This 
interpretation, first suggested by Dyson j8] for the unitary, symplectic and orthogonal 
matrix ensembles, relies on rewriting | Aat(^) | 2 as exp^X^log \zi — Zj\j- Clearly, the 
integral (12.51) looks then exactly as the partition function of the 2D Coulomb plasma 
(often called the Dyson gas) in the external field: 

Z N = J e -m*i,~*») Y[ d 2 Zj , (2.9) 

where 

Here (3 plays the role of inverse temperature (in (12. 5p f3 = 1). The first sum is the 
Coulomb interaction energy, the second one is the energy due to the external field. For 
the Hermitian and unitary ensembles the charges are confined to lines of dimension 1 
(the real line or the unit circle) but still interact as 2D Coulomb charges. So, the Dyson 
gas picture for ensembles of matrices with general complex eigenvalues distributed on the 
plane looks even more natural. It becomes especially helpful in the large N limit, where 
it allows one to apply thermodynamical arguments. 

The Dyson gas picture prompts to consider more general ensembles with arbitrary 
values of (3 ( "/3-ensembles" ) : 




(2.11) 



In general they can not be defined through matrix integrals. As we shall see, the leading 
large N contribution has a simple regular dependence on (3. However, the sub-leading 
corrections may depend on (3 in a rather non-trivial way. 

3 Exact relations at finite N 

Here we present some general exact relations for correlation functions valid for any values 
of N and (3. 
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3.1 Correlation functions: general relations 



The main objects of interest are correlation functions, i.e., mean values of functions of 
matrices. We shall consider functions that depend on eigenvalues only - for example, 
traces tr/($) = YUf( z i)- Here, /($) = /($,$ t ) is any function of $, $ f which is 
regarded as the function f(zi) = f(zi,Zi) of the complex argument Z; t (and (For the 
abuse of notation, in case of arbitrary (3 we write tr/ := Y^if{ z i) although a matrix 
realization may be not available.) Typical correlation functions which we are going to 
study are mean values of products of traces: (tr /($)), (tr /i($) tr / 2 ( < ^ ) )) and so on. 
Clearly, they are represented as integrals over eigenvalues. For instance, 

1 r / N \ N 

Z N J \i =1 ) j=1 

A particularly important example is the density function defined as 

p(z) = h^ 2 \z- Zj ), (3.1) 

3 

where S^ 2 '(z) is the 2D 5-function. Note that in our units h has dimension of [length] 2 
and p(z) is dimensionless. As it immediately follows from the definition, any correlator 
of traces is expressed through correlators of p: 

„ n 

(tr/: . . . tr/„) = h~ n j (p( Zl ) . ..p(z n )) h{z x ) . . . f n (z n ) [I d 2 Zj . (3.2) 

j'=i 

Instead of correlations of density it is often convenient to consider correlations of the field 

<p(z) = -ph'£\og\z-z j \ 2 (3.3) 

3 

from which the correlations of density can be found by means of the relation 

An(3p{z) = -A<p(z), (3.4) 

where A = Add is the Laplace operator. Clearly, tp is the 2D Coulomb potential created 
by the eigenvalues (charges). 

Handling with multi-point correlation functions, it is customary to pass to their con- 
nected parts. For example, in the case of 2-point functions, the connected correlation 
function is defined as 

(p(z 1 )p{z 2 )) c = (p{z l )p(z 2 )) - (p(zi)> (p(z 2 )) . 

The connected multi-trace correlators are expressed through the connected density cor- 
relators by the same formula (13.21) with (p(z\) . . . p{z n )) c in the r.h.s. The connected part 
of the (n + l)-point density correlation function is given by the linear response of the 
n-point one to a small variation of the potential. More precisely, the following variational 
formulas hold true: 

{p{z)) = h 1WW ' M*M*)>e = h = h 6w( Zl )6w(z 2 y (3 ' 5) 

Connected multi-point correlators are higher variational derivatives of log-Z^. These 
formulas follow from the fact that variation of the partition function over a general 
potential W inserts J2i^ 2 K z ~ z i) the integral. Basically, they are linear response 
relations used in the Coulomb gas theory 
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3.2 Loop equations 

The standard source of exact relations for correlation functions is the formal identity 

j 



E/^faW)^*)n^ = °> (3-6) 



where e({zj}) is any function of coordinates Zj bounded at infinity and E is given by 
(I2.10p . Introducing, if necessary, a cutoff at infinity one sees that the 2D integral over 
Zi can be transformed, by virtue of the Green's theorem, into a contour integral around 
infinity and so it does vanish. Being expressed in terms of correlation functions of local 
fields (such as p or (p), this identity yields, with a suitable choice of e, certain exact 
relations between them. For historical reasons, they are referred to as loop equations. 

Let us take 

««> = ™, (3.7) 

z z^ 

where X({zi}) is any symmetric function of z\, z%, . . . , zjy, then the identity reads 



E 



'-pd z .E 1 \ d g .X 
1 +t ^)X+ ' 



Z Zi [z Zip' J Z Z{ 



The singularity at the point z does not destroy the identity since its contribution is 
proportional to the vanishing integral §dzi/{z,i — z) over a small contour encircling z. 
Plugging d Zi E in the first term and using the bracket notation for the mean value, we 
have: 

h^^^ + %(z-z t )(z l -z 3 ) + ^lz~^) X + ^z^J = °- 

The second sum can be transformed by means of the simple algebraic identity 

1 



? (Z - Zi ){z - Zj ) ]^(Z- Zi){ Zi - Zj) + ^ (z - Zi f 



Finally, we arrive at the relation 



(W + E^)=0- (3-8) 

\ i Z Zi I 

where we have introduced the special notation 

or, in terms of the fields p, ip (13.1 1) . (13.31) . 

2f3h 2 L{z) = 2f3 J dW ^)p(Od 2 Z + ^ z )f + ( 2 -P)hd 2 p(z). (3.10) 

This quantity plays an important role. It is to be compared with the stress energy tensor 
in 2D CFT. 
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We call (13.81) the generating loop equation. It generates an infinite hierarchy of 
identities obeyed by correlation functions. The simplest one is obtained at X = 1: 
(L(z)) = 0. It reads: 



1 rdW(()(A<p(0) 



2nJ C-z 



d 2 ( + {{d V {z)f) + (2-0)h (&<p{z)) = 0. (3.11) 



This identity gives an exact relation between one- and two-point correlation functions 
because the mean value ((df(z)) 2 ) can be reproduced from the two-point correlation 
function (dip(z) dip(z')) by averaging over all possible directions of approaching z' — > z: 

({d^{z)f) = Jim i- ^ (dp(z) d<p{z+ee i6 )) d9 . (3.12) 

Another interesting choice is X = <p((), where ( ^ z, then d Zi X = f3K/(£ — Zj) and (13.81) 
yields 

(L(.-M0 + ^:^ (C) )=0. (3.13) 

Acting by d^, we get (^L(z)Bip(()' K J = ir(3(p(()} /(( — z). Further, acting by to both 
sides, we obtain the relation 

(iW/)(c) > = w^l> ( 3. 14 ) 
{z-cy z-( 

A similar but longer calculation gives a generalization of this identity for m-point func- 
tions: 

m ( 1 1 \ 

(L(z)p(d) . . . p(c m )) = E r^7V2 + —T7 ^ MCO • • • • ( 3 - 15 ) 

i=l sij z Si ) 

It has the form of the conformal Ward identity for primary field of conformal dimensions 
1, with L(z) playing the role of the (holomorphic component of) the stress energy tensor 
in CFT. 



4 Large N limit 

Starting from this section, we study the large A" limit 

N^oo, h^O, %N = t finite. (4.1) 

We shall see that in this limit meaningful analytic and algebro-geometric structures 
emerge, as well as important applications in physics. 



4.1 Solution to the loop equation in the leading order 

It is instructive to think about the large A" limit under consideration in terms of the Dyson 
gas picture. Then the limit we are interested in corresponds to a very low temperature 
of the gas, when fluctuations around equilibrium positions of the charges are negligible. 
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The main contribution to the partition function then comes from a configuration, where 
the charges are "frozen" at their equilibrium positions. It is also important that the 
temperature tends to zero simultaneously with increasing the number of charges, so the 
plasma can be regarded as a continuous fluid at static equilibrium. Mathematically, all 
this means that the integral is evaluated by the saddle point method, with only the 
leading contribution being taken into account. As H — > 0, correlation functions take their 
"classical" values ((p(z)) = <fo{z), (p(z)) = Po{z), and multi-point functions factorize in 
the leading order: (dcp(z) d(p(z')) = d(po(z)d(po(z'), etc. Then the loop equation (13.111) 
becomes a closed relation for <p : 



1 r dW(C)Aip (0 J2/ 



27T J C 



d\ + (dMz)) = 0, (4.2) 



where we have ignored the last term because it has a higher order in h. Applying 8 to the 
both terms, we get: — dW(z)A(p (z) + d<p (z)A(p (z) = 0. Since A<p (z) = — 4Tr/3p (z) 
(see (13.41) ). we obtain 

p (z) [d<p Q (z)-dW{z)]=0. (4.3) 

This equation should be solved with the additional constraints / po(z)d 2 z = t (normal- 
ization) and po(z) > 0. The equation tells us that either d(po(z) = dW(z) or po(z) = 0. 
Applying d to the former, we get Aip (z) = AW(z). This gives the solution for p : 

p J z ) = - AW ( Z ) in the bulk (4.4) 

Here, "in the bulk" means "in the region where po > 0". The physical meaning of 
the equation d<p (z) = dW(z) is clear. It is just the condition that the charges are in 
equilibrium (the saddle point for the integral). Indeed, the equation states that the total 
force experienced by a charge at any point z where po 7^ is zero, i.e., the interaction 
with the other charges, dipo(z), is compensated by the force dW(z) due to the external 
field. 



4.2 Support of eigenvalues 

Let us assume that 

a(z) := --AW(z) > 0. (4.5) 

47T 

For quasiharmonic potentials, a(z) = If, according to (I4.4p . p = <r//3 everywhere, 
the normalization condition for p in general can not be satisfied. So we conclude that 
Po = o/(3 in a compact bounded domain (or domains) only, and outside this domain one 
should switch to the other solution of (14.31) . po = 0. The domain D where po > is called 
support of eigenvalues or droplet of eigenvalues. In general, it may consist of several 
disconnected components. The complement to the support of eigenvalues, D c = C \ D, is 
an unbounded domain in the complex plane. For quasiharmonic potentials, the result is 
especially simple: po is constant in D and in D c . 

To find the shape of D is a much more serious problem. It appears to be equivalent 
to the inverse potential problem in 2D. The shape of D is determined by the condition 
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Figure 1: The normal displacement of the boundary. 



d(fio{z) — dW(z) (valid at all points z inside D) and by the normalization condition. One 
can write them in the form 



47T Jd z — ( 



a(()d 2 ( = f3t . 



dW{z) 



for all z G D 



(4.6) 



The integral over D in the first equation can be transformed to a contour integral by 
means of the Cauchy formula. As a result, the first equation reads: 

r dW(Qdt 

TdD Z — ( 







for all z G D. 



(4.7) 



This means that the domain D has the following property: the function dW(z) on its 
boundary is the boundary value of an analytic function in its complement D c . 

The connection with the inverse potential problem is most straightforward in the 
quasiharmonic case, where dW(z) = —z + V'(z), AW(z) = AddW(z) = —4. The nor- 
malization then means that the area of D is equal to /3irt. Assume that: i) V(z) = J2k tkZ k 
is regular in D (say a polynomial), ii) e D (it is always the case when W has a maximum 

at 0), iii) D is connected. Then equation (14.71) acquires the form / 5 ~ = V'( z ) 

2m J do C — z 

for zGD. Expanding it near z = 0, we get: 



CC k dC = - 



C k d 2 ( 



(4- 



2nik JdD 3 " nk Jd 

We see that the "coupling constants" tk are harmonic moments of D c = C \ D and the 
area of D is TT/3t. It is the subject of the inverse potential problem to reconstruct the 
domain from its area and harmonic moments. In general, the solution is not unique. But 
it is known that locally, i.e., for a small enough change t —>■ t + 5t, t& — > + Stk there is 
only one solution. 

As an explicitly solvable example, consider the Ginibre-Girko ensemble with the parti- 
tion function (12. 8p . In this case the support of eigenvalues is an ellipse with the half-axes 



t 



l + 2|t 2 



\ 1 — 2I* 2 



t 



2 It, 



\ 1 + 2|* 2 
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centered at the point z G = — 1 2 * , and with the angle between the big axis and the 

l-4|t 2 | 2 

real line being equal to ^argt 2 - For the model with the potential W = —\z\ 2 + 23^(t 3 ^ 3 ) 
the droplet of eigenvalues is bounded by a hypotrochoid [i~2]. 



4.3 Small deformations of the support of eigenvalues 

Coming back to the general case, let us examine how the shape of D changes under a 
small change of the potential W with t = hN fixed. It is convenient to describe small 
deformations D — ► D, by the normal displacement 5n(() at a boundary point £ (Fig. [p. 
Consider a small variation of the potential W in the condition (14. TJ) . To take into account 
the deformation of the domain, we write, for any fixed function /, 
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<f ftt)d() = <f /(CR = 2i/ df(()d 2 ( « 2z/ 5/(C)KC)|dCI 

/ Jd(SD) J 8D JdD 



ld(8D) 

(here 5D = D \ D) and thus obtain from (14.7ft : 



9D -2 — C 2 7aD -2 — C 

This integral equation for 5n(Q can be solved in terms of the exterior Dirichlet boundary 
value problem. Given any smooth function f(z), let f H (z) be its harmonic continuation 
from the boundary of D to its exterior, i.e., a unique function such that Af H = in D c 
and regular at oo, and f H {z) = f(z) for all z G dD. Explicitly, a harmonic function can 
be reconstructed from its boundary value by means of the formula 

f H (z) = -^{ f(0d n G(z,0m. (4.10) 

Z7T J dD 

(Here and below, d n is the normal derivative at the boundary, with the outward pointing 
normal vector.) The main ingredient of this formula is the Green's function G(z,£) of the 
domain D c characterized by the properties A z G(z, () = 2nS^(z — Q in D c , G(z,Q = 
if z or ( G dD. As ( — > z, it has the logarithmic singularity G(z, () — > log \z — C\. 

r d(5W H (0)d£ 

Consider the integral i> which is obviously equal to for all z inside 

JdD Z — 4 

D, subtract it from the first term in (14. 9 p and rewrite the latter as an integral over the 
line element \d£\. After this transformation (14.91) acquires the form 

— / \d£\ = for all z G D, (4.11) 

2ni JdD z — 4 

where 

R(z) = AW(z)5n(z) + d~ (8W(z) - (5W) H (z)) (4.12) 

is a real- valued function on the boundary contour <9D. The superscript indicates that the 
derivative is taken in the exterior of the boundary. By properties of Cauchy integrals, 
it follows from (14.111) that i?(£)/r(£), where r(£) = d£/\d£\ is the unit tangential vector 
to the boundary curve, is the boundary value of an analytic function h(z) in D c such 
that h(oo) = 0. For z G D c , this function is just given by the integral in the l.h.s. of 
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(14. lip . Variation of the normalization condition (the second equation in (14. 6p ) yields, in 
a similar manner: 

/ R(£M\ = (4.13) 

JdD 

This relation implies that the zero at oo is at least of the 2-nd order. 

The following simple argument shows that an analytic function with these properties 
must be identically zero. Let w(z) be the conformal map from D c onto the unit disk such 
that tu(oo) = and the derivative at oo is real positive. By the well known property of 
conformal maps we have 

d z e iaxgw'(z) _ dw 



\dz\ \dw\ 

along the boundary curve. Therefore, r(z) = i\w'(z) \ w(z)/w'(z) and we thus see that 

R(z) w'(z) 
i\w'(z)\w(z) 

is the boundary value of the holomorphic function h(z). Since w'(z) ^ in D c , the 
function g(z) = h(z) w(z)/w'(z) is holomorphic there with the purely imaginary boundary 
value ■ R ^\, . Then the real part of this function is harmonic and bounded in D c and is 
equal to on the boundary. By uniqueness of a solution to the Dirichlet boundary value 
problem, dtg(z) must be equal to identically. Therefore, g(z) takes purely imaginary 
values everywhere in D c and so is a constant. By virtute of condition (14.131) this constant 
must be which means that R(z) = 0. 

Therefore, one obtains the following result for the normal displacement of the bound- 
ary caused by a small change of the potential W — > W + 5W: 

(4,4, 



4.4 From the support of eigenvalues to an algebraic curve 

There is an interesting algebraic geometry behind the large N limit of matrix models. 
For simplicity, here we consider models with quasiharmonic potentials. 

The boundary of the support of eigenvalues is a closed curve in the plane without 
self- intersections. If V'(z) is a rational function, then this curve is a real section of a 
complex algebraic curve of finite genus. In fact, this curve encodes the 1/N expansion 
of the model. In the context of Hermitian 2-matrix model such a curve was introduced 
and studied in [231124]. 

To explain how the curve comes into play, we start from the equation d(po = dW, 

1 r d 2 ( 

which can be written in the form z — Viz) = G(z) for zgD, where G(z) = — -. 

7T Jd z — Q 

Clearly, this function is analytic in D c . At the same time, V'(z) is analytic in D and 
all its singularities in D c are poles. Set S(z) = V'(z) + G(z). Then S(z) = z on the 
boundary of the support of eigenvalues. So, S(z) is the analytic continuation of z away 
from the boundary. Assuming that poles of V are not too close to <9D, S(z) is well- 
defined at least in a piece of D c adjacent to the boundary. The complex conjugation 
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yields S(z) = z, so the function S(z) = S(z) must be inverse to the S(z): S(S(z)) = z 
("unitarity condition"). The function S(z) is called the Schwarz function [25] . 

Under our assumptions, S(z) is an algebraic function, i.e., it obeys a polynomial 
equation P(z, S(z)) = of the form 

d+l 

P{z,S(z))= a nl z n (S(z)) l = 0, 

n,l=l 

where aj^ = a n \ and d is the number of poles of V'(z) (counted with their multiplicities). 
Here is a sketch of proof. Consider the Riemann surface £ = D c U <9D U (D c )* (the 
Schottky double of D c ). Here, (D c )* is another copy of D c , with the local coordinate z, 
attached to it along the boundary. On A7, there exists an anti-holomorphic involution 
that interchanges the two copies of D c leaving the points of dD fixed. The functions z 
and S(z) are analytically extendable to (D c )* as S(z) and z respectively. We have two 
meromorphic functions, each with d+l poles, on a closed Riemann surface. Therefore, 
they are connected by a polynomial equation of degree d + 1 in each variable. Hermiticity 
of the coefficients follows from the unitarity condition. 

The polynomial equation P(z, z) = defines a complex curve r with anti-holomorphic 
involution (z, z) i— >■ (z, z). The real section is the set of points such that z — z. It is the 
boundary of the support of eigenvalues. 

It is important to note that for models with non-Gaussian weights (in particular, with 
polynomial potentials of degree greater than two) the curve has a number of singular 
points, although the Riemann surface £ (the Schottky double) is smooth. Generically, 
these are double points, i.e., the points where the curve crosses itself. In our case, a 
double point is a point 6 D c such that S(z^) = z^ but z^ does not belong to the 
boundary of D. Indeed, this condition means that two different points of S, connected 
by the antiholomorphic involution, are stuck together on the curve r, which means 
self-intersection. The double points play the key role in deriving the nonperturbative 
(instanton) corrections to the large N matrix models results (see [26J for details). 

4.5 Free energy and correlation functions 

Free energy. The leading contribution to the free energy in the large N limit is 

F = lim (h 2 log Z N ) . 

It is determined by the maximal value of the integrand in (12. 9p . i.e., by the extremum of 
the function E. In the continuous approximation, one can represent it as a functional of 
the density: 

- f3h 2 E[p] =(3 J J p(z)p(() log \z - (\d 2 zd 2 ( + J W(z)p(z)d 2 z (4.15) 

and find its minimum with the constraint / p d 2 z = t. Introducing a Lagrange multiplier 
A, we get the equation 

2(3 J log \z - (\p(()d 2 ( + W{z) + A = 
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which is solved, as expected, by the function p discussed in section 4.1. Assuming that 
W(0) = 0, we find A = <p (0) and 



1 



a(z) lot 



1 1 



a{Qd 2 zd\. 



(4.16) 



Some results for correlation functions. Here are some results for the correlation 
functions obtained by the variational technique. (For details of the derivation see |lUj.) 
They are correct at distances much larger than the mean distance between the charges. 
The leading contribution to the one-trace function was already found in Section 4.1: 



(4.17) 



(tr/(d>)) = —j D a(z)f(z)d 2 z + 0(l). 



The connected two-trace function is: 

(tr /($) tr g(d>)) c = -L J VfVgd 2 z - -L f fd n g H \dz\ + 0(h) (4.18) 



(here V/Vp = d x fd x g + d y fd y g). In particular, for the connected correlation functions 
of the fields <p(z) this formula gives: 



<<p(*M*0>c 



-2ptr log 



+ Oih 6 ) inside D 



2(3h 2 G(z, z')-G(z, oo)-G(z', oo) -log 



\z — z 



-0(h 3 ) outside D 



where G is the Green's function of the Dirichlet boundary value problem and 



exp 



lim(log|£| + G(£,oo)) 



(4.19) 



(4.20) 



is the external conformal radius of the domain D (the Robin's constant). This result 
is valid if \z — z'\ ^> V%. The 2-trace functions are universal, i.e., they depend on the 
shape of the support of eigenvalues only and do not depend on the potential W explicitly. 
They resemble the two-point functions of the Hermitian 2-matrix model found in [27J; 
they were also obtained in [28] in the study of thermal fluctuations of a confined 2D 
Coulomb gas. The structure of the formulas indicates that there are local correlations in 
the bulk as well as strong long range correlations at the edge of the support of eigenvalues. 
(See [29J for a similar result in the context of classical Coulomb systems). 



5 The matrix model as a growth problem 
5.1 Growth of the support of eigenvalues 

When N increases at a fixed potential W, one may say that the support of eigenvalues 
grows. More precisely, we are going to find how the shape of the support of eigenvalues 
changes under t — > t + St, where t = Nh, if W stays fixed. 
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The starting point is the same as for the variations of the potential, and the calcula- 
tions are very similar as well. Variation of the conditions (14.61) yields 



JdD 



AW(QSn(Q 
z-C 



\d(\=0 for all z G D, / AW{OSn{0 \d(\ = -AnpSt . 

JdD 



The first equation means that AW(z)Sn(z)^jj- is the boundary value of an analytic 
function h(z) such that h(z) = —47rf3St/z + 0(z~ 2 ) as z —>■ oo. The solution for the 
Sn(z) is again expressed in terms of the Green's function in D c : 



Sn(z) 



pst 

l-naiz 



d n G(oo, z). 



(5.1) 



For quasiharmonic potentials (with o = 1/vr), the formula simplifies: 

Sn(z) = — — St d n G(oo, z). 



(5.2) 



Identifying t with time, one can say that the normal velocity of the boundary at any 
point z, V n (z) = Sn(z)/8t, is proportional to gradient of the Green's function: V n (z) oc 
-d n G(oo,z). 

If the domain is connected, G can be expressed through the conformal map w(z) from 
D c onto the exterior of the unit circle: 



G(z 1 ,z 2 ) = log 



w(zi) - w(z 2 ) 



w(z 1 )w(z 2 ) 



(5.3) 



In particular, G(oo,z) 



\og\w(z)\. As \z\ — > oo, w(z) = z/r + 0(1), where r is the 



Robin's constant which enters eq. (14. 19j) . It is easy to see that <9„log \ w(z) 
dD, so one can rewrite the growth law (I5.2p as 5n(z) = | |u/( z )l St. 



\w \ z)\ on 



5.2 Laplacian growth 

The growth law (15.21) is common to many important problems in physics. The class 
of growth processes, in which dynamics of a moving front (an interface) between two 
distinct phases is driven by a harmonic scalar field is known under the name Laplacian 
growth. The most known examples are viscous flows in the Hele-Shaw cell, filtration 
processes in porous media, electrodeposition and solidification of undercooled liquids. A 
comprehensive list of relevant papers published prior to 1998 can be found in [IB] . 

Viscous flows in Hele-Shaw cell. To be specific, we shall speak about an interface 
between two incompressible fluids with very different viscosities on the plane (say, oil 
and water). In practice, the 2D geometry is realized in the Hele-Shaw cell - a narrow 
gap between two parallel glass plates. For a review, see [TH [15] . The velocity field in a 
viscous fluid in the Hele-Shaw cell is proportional to the gradient of pressure p (Darcy's 
law): 

b 2 

V = -KVp, K = —. 
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water 




Figure 2: The Hele-Shaw cell. 



The constant K is called the filtration coefficient, /i is viscosity and b is the size of the gap 
between the two plates. Note that if jj, — > 0, then Vp — > 0, i.e., pressure in a fluid with 
negligibly small viscosity is uniform. Incompressibility of the fluids (VV = 0) implies 
that the pressure field is harmonic: Ap = 0. By continuity, the velocity of the interface 
between the two fluids is proportional to the normal derivative of the pressure field on 
the boundary: V n = —Kd n p. 

To be definite, we assume that the Hele-Shaw cell contains a bounded droplet of 
water surrounded by an infinite "sea" of oil (another possible experimental set-up is an 
air bubble surrounded by oil or water). Water is injected into the droplet while oil is 
withdrawn at infinity at a constant rate, as is shown schematically in Fig. [21 The latter 
means that the pressure field behaves as p oc — log \z\ at large distances. We also assume 
that the interface is a smooth closed curve T. As it is mentioned above, one may set p = 
inside the water droplet. However, pressure usually has a jump across the interface, so p 
in general does not tend to zero if one approaches the boundary from outside. This effect 
is due to surface tension. It is hard to give realistic estimates of the surface tension effect 
from fundamental principles, so one often employs certain ad hoc assumptions. The most 
popular one is to say that the pressure jump is proportional to the local curvature, k, of 
the interface. 

To summarize, the mathematical setting of the Saffman- Taylor problem is as follows: 



V n = -d n p 
Ap = 
p — > — log \z\ 
p = 

— <pi ) — — 



UK 



on r 

in oil 

in oil as \z\ 
in water 
across T 



■DO 



(5.4) 



Here v is the surface tension coefficient. (The filtration coefficient is set to be 1.) The 
experimental evidence suggests that when v is small enough, the dynamics becomes 
unstable. Any initial domain develops an unstable fingering pattern. The fingers split into 
new ones, and after a long lapse of time the water droplet attains a fractal-like structure. 
This phenomenon is similar to the formation of fractal patterns in the diffusion-limited 
aggregation. 

Comparing (15.21) and (15.41) . we identify D and D c with the domains occupied by 
water and oil respectively, and conclude that the growth laws are identical, with the 
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pressure field being given by the Green function: p(z) = G(oa,z), and p = on the 
interface. The latter means that supports of eigenvalues grow according to (15.41) with 
zero surface tension, i.e., with v = in (15 .4p . Neglecting the surface tension effects, one 
obtains a good approximation unless the curvature of the interface becomes large. We 
see that the idealized Laplacian growth problem, i.e., the one with zero surface tension, 
is mathematically equivalent to the growth of the support of eigenvalues in ensembles of 
random matrices M or C and in general /^-ensembles. 

Finite-time singularities as critical points. As a matter of fact, the Laplacian 
growth problem with zero surface tension is ill-posed since an initially smooth interface 
often becomes singular in the process of evolution, and the solution blows up. The role of 
surface tension is to inhibit a limitless increase of the interface curvature. In the absence 
of such a cutoff, the tip of the most rapidly growing finger typically grows to a singularity 
(a cusp). In particular, a singularity necessarily occurs for any initial interface that is 
the image of the unit circle under a rational conformal map, with the only exception of 
an ellipse. 

An important fact is that the cusp-like singularity occurs at a finite time t = t c , i.e., 
at a finite area of the droplet. It can be shown that the conformal radius of the droplet 
r (as well as some other geometric parameters), as t — > t c , exhibits a singular behavior 

r - r c oc (t c - t)~ 7 

characterized by a critical exponent 7. The generic singularity is the cusp (2,3), which 
in suitable local coordinates looks like y 2 = x 3 . In this case 7 = — ~. 

A similar phenomenon has been known in the theory of random matrices for quite a 
long time, and in fact it was the key to their applications to 2D quantum gravity and 
string theory. In the large N limit, the random matrix models have critical points - the 
points where the free energy is not analytic as a function of a coupling constant. As we 
have seen, the Laplacian growth time t should be identified with a coupling constant of 
the normal or complex matrix model. In a vicinity of a critical point, 

F ~F rcg +a (t c -t) 2 -\ 

where the critical index 7 (often denoted by 7 str in applications to string theory) depends 
on the type of the critical point. Accordingly, the singularities show up in correlation 
functions. Using the equivalence established above, we can say that the finite-time blow- 
up (a cusp-like singularity) of the Laplacian growth with zero surface tension is a critical 
point of the normal and complex matrix models. Remarkably, the evolution can be 
continued to the post-critical regime as a dynamics of "shock lines" [30J . 
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